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Abstract 

We present results from the simulation of a two-coupling spin-1 model 
with states 0, ±1 and nearest neighbour interaction. By a suitable 
choice of couplings we are able to drastically reduce the effects of 
corrections to scaling. Our estimates for the critical exponents are 
v = 0.6299(3) and 77 = 0.0359(10). For the universal ratio Q = 
(m 2 ) 2 /(m 4 ) we obtain Q = 0.6240(2). The universal ratio of par- 
tition functions with antiperio die/periodic boundary conditions, re- 
spectively, is Z a /Z p = 0.5425(2). 



1 Introduction 



The determination of parameters like critical exponents from finite size scal- 
ing is plagued by the appearance of corrections to scaling. Consider, e.g., 
an Ising model on a cubic lattice of linear extension L. The Binder cumulant 

(m 4 ) 

behaves at the critical point like 

Q(L) = Q* + qL- UJ + ... (2) 

m denotes the magnetization per spin. Q* is universal, u ~ 0.8 denotes the 
leading correction to scaling exponent. Many more terms appear in eq. (0), 
governed by subleading exponents and combinations of them. If they have 
non-negligible coefficients in front of them, they can hamper or even make 
impossible a reliable determination of the infinite L behaviour. 

In the language of the renormalization group the presence of strong cor- 
rections to scaling is related to the fact that one is far away from the renor- 
malization group fixed point. 

In this article, we summarize the results of our attempt to improve the 
scaling properties of the 3D Ising model by suitably tuning the parameters 
in a two-coupling spin-1 model that has already been studied by Blote et 
al. 0. 

We will be mainly interested in the determination of the exponents v and 
rj that will be obtained from the scaling laws (i = 1,2) 

|| = a i L 1 /*(l + & i L- + ...) (3) 

and 

X = c + L 2 -v(l + dL- UJ + ...). (4) 

Here, R2 denotes the cumulant introduced above, and Ri is the ratio of par- 
tition functions Z a /Z p , where the subscripts label periodic and antiperiodic 
boundary conditions, respectively. Antiperiodic boundary conditions are ap- 
plied in one of the three lattice directions only. The susceptibility \ is defined 
by 

X = L 3 (m 2 ) . (5) 
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2 Improved Spin-1 Ising Model 



Following |2j, we consider a spin-1 Ising model on a simple cubic lattice 
defined through 

S = -P £ s lS , + Dj2sl (6) 

<i,j> i 

The Sj take values 0, ±1, and the spin-spin interaction is a sum over all near- 
est neighbour pairs. The Boltzmann factor is given by exp(— S). Blote et 
al. chose a fixed D = In 2, and varied (3 to tune to criticality. We followed a 
different procedure, to be described in some detail elsewhere. From prior sim- 
ulations ||, estimates of the universal quantities R* were available, namely 
.R* = 0.5425(10) and R* 2 = 0.624(1). For a number of lattice sizes Lj rang- 
ing from 3 to 14 a sequence of couplings (Pj, Dj) was determined such that 
Ri(j3j, Dj, Lj) = R* within the statistical accuracy. As L increases the pair 
(Pj, Dj) converges to a point where leading order corrections to scaling vanish. 
From our numerical results we extrapolate to (/3oo,Doo) = (0.3832,0.6241). 
We found that this point is approached on the line 

D ~ —2.04 + 6.95 (3 . (7) 

It turns out that eq. (0) is an approximation to the critical line. It is reas- 
suring that the estimate p c = 0.3934214(8), together with D = In 2 in 
is close to this line. Our next goal was to approximately identify a second 
one-dimensional submanyfold of the 2-coupling space, namely that where 
the leading correction to scaling vanishes. This was possible by fitting the 
derivatives of the Ri(f3j, Dj, Lj) at (Poo, A») with respect to the two coupling 
parameters to a certain linear equation. Details will be reported in JJ. Our 
fit results suggested that the critical line should be approached by varying P 
while adjusting D according to 

D = 3 (p - 0.3832) + 0.6241 . (8) 

A first estimate for p c was then obtained applying the well established method 
looking for the crossings of Ri(L) with Ri(2L). We used lattices of size 
L = 4,8,16,32 and obtained p c = 0.383245(10), with D given by eq. (g), 
i.e. D c = 0.624235. This is significantly different from the D = In 2 used in 
ref. §. 
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3 Simulation Results 



Monte Carlo simulations were now performed at (3 C = 0.383245, fixing D 
according to eq. (H). We simulated on cubic lattices with linear extension 
L ranging from 4 to 56. State of the art cluster algorithms were employed. 
Two measurements were separated by three single cluster updates. For the 
computation of Z a /Z p a variant of the boundary flip cluster algorithm Q 
was implemented. The runs for the determination of the critical parameters 
described in the previous section took about 3 months of CPU on Pentium 
166-MMX PCs, while the final production runs consumed a total of 1 year. 

The data relevant for the present article are summarized in tables |Ij and 
Table |l] gives the quantities Ri, i — 1, 2 and xjL 2 . In addition, we quote the 
number of measurements in the last column. In table |2] we give the estimates 
of the partial and total derivatives of the Ri with respect to (3. Note that 
eq. (H) implies 

dR t _ dRj 8Rj 
dp d/3 3D' 

The quantities in table |2] are multiplied by a factor f(L) = L -1 / - 6 ^ which 
compensates to good precision for the leading divergent behaviour with in- 
creasing lattice size. 

By a reweighting technique we have access also to the same set of ob- 
servables at a set of four neighbouring P values, ranging from 0.383225 to 
0.383265 in steps of 0.000010. Having control over a neighbourhood of the 
simulation point is important, since we started without a high precision es- 
timate of (3 C . 



4 Fitting the Data 
4.1 Fitting Z a /Z p , Q 

We first fitted the Ri in order to obtain estimates for the universal quantities 
R*. The scaling law in question is eq. (^). However, looking carefully at 
table [I] reveals that corrections to scaling are small. We therefore fitted the 
data with the law 

Ri(L, p MC ) = R* + -j7f{L, Pmc) A/5 . (9) 
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L 


Z IZ 


O 


y/L 2 

A/ ^ 


stat/10 6 


4 


.53599(07) 


.62408(05) 


.87609(13) 


30 


6 


.54080(07) 


.62491(06) 


.88104(14) 


60 


8 


.54193(09) 


.62460(07) 


.87706(16) 


36 


10 


.54209(18) 


.62464(13) 


.87272(33) 


10 


12 


.54248(19) 


.62416(14) 


.86778(35) 


10 


14 


.54230(21) 


.62426(15) 


.86399(38) 


10 


16 


.54263(17) 


.62409(13) 


.86000(32) 


13 


18 


.54251(23) 


.62410(17) 


.85697(42) 


8 


20 


.54286(28) 


.62379(21) 


.85332(50) 


3 


22 


.54244(30) 


.62404(23) 


.85115(55) 


7 


24 


.54229(23) 


.62402(17) 


.84914(43) 


3 


28 


.54245(31) 


.62386(23) 


.84373(55) 


6 


32 


.54259(31) 


.62390(23) 


.84010(56) 


6 


36 


.54213(41) 


.62420(30) 


.83742(73) 


5 


40 


.54227(32) 


.62403(23) 


.83392(56) 


7 


48 


.54250(39) 


.62395(29) 


.82803(69) 


5 


56 


.54272(63) 


.62382(47) 


.82351(12) 


2 



Table 1: Results for i?i = Z a /Z p , R2 = Q and the susceptibility- 
divided by L 2 at = 0.383245. The last column gives the number of 
measurements in units of 10 6 . 
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L 


f(L) dRx/dp 


f(L) dRx/dfi 


f(L) 8R 2 /d(3 


f(L) dR 2 /d(3 


4 


-1.09974(17) 


-.63210(11) 


.65894(15) 


-.09882(03) 


6 


-1.13595(17) 


-.64829(10) 


.66236(15) 


-.09685(03) 


8 


-1.14621(20) 


-.65290(12) 


.66295(18) 


-.09614(03) 


10 


-1.15051(41) 


-.65498(24) 


.66276(37) 


-.09577(07) 


12 


-1.15293(45) 


-.65567(25) 


.66307(40) 


-.09571(07) 


14 


-1.15385(49) 


-.65661(27) 


.66310(46) 


-.09555(08) 


16 


-1.15432(44) 


-.65652(25) 


.66290(39) 


-.09554(07) 


18 


-1.15547(59) 


-.65692(33) 


.66294(55) 


-.09545(11) 


20 


-1.15603(71) 


-.65732(39) 


.66322(64) 


-.09553(11) 


22 


-1.15542(79) 


-.65714(43) 


.66310(74) 


-.09574(17) 


24 


-1.15675(61) 


-.65790(34) 


.66419(55) 


-.09568(11) 


28 


-1.15565(80) 


-.65730(44) 


.66199(75) 


-.09520(20) 


32 


-1.15664(82) 


-.65793(46) 


.66283(75) 


-.09523(18) 


36 


-1.1568(12) 


-.65776(61) 


.6638(11) 


-.09525(49) 


40 


-1.15665(84) 


-.65785(47) 


.66326(84) 


-.09557(37) 


48 


-1.1563(11) 


-.65752(59) 


.6635(11) 


-.09513(58) 


56 


-1.1562(17) 


-.65739(97) 


.6628(19) 


-.0958(12) 



Table 2: Results for partial and total derivatives of R\ = Z a /Z 1 
R 2 = Q, multiplied by the factor f(L) = L" 1 / - 63 , at (3 = 0.383245. 
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-^min 


R[ 




X 2 /dof 


6 


0.54143(5) 


0.3832516(8) 


15.49 


8 


0.54213(8) 


0.3832470(8) 


1.78 


10 


0.54240(10) 


0.3832453(9) 


0.79 


12 


0.54251(11) 


0.3832447(9) 


0.49 


16 


0.54260(15) 


0.3832442(11) 


0.46 


20 


0.54252(25) 


0.3832446(14) 


0.55 


24 


0.54226(37) 


0.3832457(20) 


0.27 


Table 3: Fitting R x with eq. © 


-^min 




& 


X 2 /clof 


6 


0.62468(4) 


0.3832523(10) 


4.59 


8 


0.62447(6) 


0.3832499(10) 


2.22 


10 


0.62429(7) 


0.3832479(12) 


1.39 


12 


0.62414(8) 


0.3832465(11) 


0.42 


16 


0.62405(12) 


0.3832457(14) 


0.31 


20 


0.62393(18) 


0.3832447(18) 


0.28 


24 


0.62400(27) 


0.3832452(24) 


0.23 



Table 4: Fitting R 2 with eq. © 



Here we have included a term which (to first order) corrects for deviations 
from being at criticality. (3mc is our simulation coupling 0.383245, and the 
dRi/d/3 are taken from table |[ The fit parameters are R* and [3 C entering 
through A/3 = /3mc — fie- To check for effects from corrections to scaling, the 
fits were done on a sequence of data sets obtained by discarding data with 
L < L min . We first fitted separately R± and R 2 . The results are reported in 
tables [3] and[|. Then we fitted all data together with three parameters (R^, 
i?2) an d Pc)- The results are presented in table ||. 

Our fit results for Q* and the /3 c 's from the Q-fits (as function of I/ min ) are 
shown in figures |I| and 0. The plots justify our final estimates R 2 = 0.6240(2) 
and f3 c = 0.383245(1). For R\ we find 0.5425(2). 
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-^min 




R\ 




X 2 /dof 


6 


0.62467(4) 


0.54142(4) 


0.3832519(5) 


10.05 


8 


0.62441(5) 


0.54206(7) 


0.3832481(6) 


2.19 


10 


0.62421(6) 


0.54231(9) 


0.3832463(7) 


1.21 


12 


0.62408(7) 


0.54245(11) 


0.3832453(8) 


0.51 


16 


0.62399(8) 


0.54254(14) 


0.3832448(9) 


0.42 


20 


0.62393(13) 


0.54252(22) 


0.3832446(11) 


0.41 


24 


0.62403(13) 


0.54230(22) 


0.3832455(11) 


0.25 



Table 5: Fitting simultaneously R<i and R\ with eq. (|S|) 



0.6248 




10 



15 



20 



25 



30 



35 



Figure 1: Fit results for R\ as function of smallest lattice in- 
cluded in the fit. 
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ft 

0.383254 | . . 1 . . 1 

0.383252 -{ 
0.383250 - 

0.383248 -I - 
0.383246 - j. . T . T ...... .....T ....... 

i-4-4-i-^--4-— ■ 

0.383244 - 1 I 

0.383242 - 

0.383240 I 1 1 1 1 1 

5 10 15 20 25 30 35 

^min 

Figure 2: Fit results for (3 C as function of smallest lattice included 
in the fit. 
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4.2 Fitting the Derivatives of Ri 



We fitted our results for the derivatives of the Ri with respect to (3, according 
to the law 

fW- do) 

Our results for 02, i.e., for the derivative of the cumulant, are presented in 
table ||. The ^-results of this table, together with the corresponding results 
for Z a / Z p are plotted in figure [| Obviously, the derivatives of the cumulant 
scale much better than those of Z a /Z p . 

For the derivative of Ri, including the leading correction to scaling plus a 
term ~ L~ x with x of order 2 improved the fits. Note, however, that there is 
more than one exponent or combination of exponents of order 2 (e.g,. 1 / v+uj 
or a/), and fits of this type do not have a sound theoretical basis. The nice 
scaling behaviour of our data (especially of the cumulant related quantities) 
allows us to avoid getting into the mess of multi-parameter fits. 

To check for the systematic dependence on the location of /3 C we repeated 
the fit for the Q-derivative on data from five shifted (3- values from 0.383243 to 
0.383247 in steps of 0.000001, covering thus two standard deviations around 
our f3 c estimate. Figure |] shows that the variation of the v estimate to this 
(3 allows us to put as our final estimate 

v = 0.6299(3). (11) 



4.3 Fitting the Susceptibility 

We fitted the susceptibility data given in table [TJ, and the corresponding 
data at shifted (3- values with eq. (f|). It turned out that the first correction 
to scaling is negligible. However, the fit results depend to some extent on 
the location of f3 c . The results for the exponent 77 for (3 C = 0.383245 are 
given in table |7|. Figure |5| shows the rj estimates from table [7| and also those 
from /5-values chosen two standard deviations above (upper data) and below 
(lower data) our (3 C estimate. 

We then tried to do the fits in an alternative way. Define a function (3 C (L) 
by requiring that for any L the relation Q*(L, (3 C (L)) = Q* = 0.6240. The 
susceptibility as a function of L and (3 C {L) should behave to leading order as 

x(L,f3 c (L))=5 + dL 2 -^. (12) 
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&% 


V 


v 2 /dof 

A, / v -*- v -'■ , ■ 


6 


0.66160(42) 


0.62969(11) 


0.79 


8 


0.66247(63) 


0.62987(14) 


0.57 


10 


0.6622(11) 


0.62982(22) 


0.60 


12 


0.6626(12) 


0.62989(24) 


0.64 


14 


0.6624(18) 


0.62986(36) 


0.70 


16 


0.6622(19) 


0.62983(37) 


0.77 


18 


0.6631(23) 


0.62998(42) 


0.84 


20 


0.6643(31) 


0.63018(56) 


0.90 


24 


0.6662(51) 


0.63051(87) 


1.18 



Table 6: Fitting the derivative of R2 with respect to (3. The fit 
function is given in eq. fllCf ). 



v 

0.632 1 1 1 1 1 1 

0.631 t 

°- 630 i f T ] \ [- V. I ' 

0.629 - i * * 

0.628 - 
0.627 
0.626 - 

0.625 1 1 1 1 1 1 

5 10 15 20 25 

^min 



Figure 3: Fit results for v from fitting dRi/d[3 with eq. (|IOl ). 
The data points with the better scaling behaviour belong to the 
cumulant R2. 
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0.6315 
0.6310 
0.6305 
0.6300 
0.6295 



fe^ ' a J.: if 




0.6290 




10 



15 
L 



20 



25 



Figure 4: Fit results for v from fitting dR 2 /d/3 with eq. (PH). 
To check for effects of errors in the location of (3 C , the //-estimate 
is given for the five /3-values 0.383243 to 0.383247 in steps of 
0.000001 (left to right). 
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Lmin 


c 


d 


V 


X 2 /dof 


4 


-0.5216(56) 


0.95751(59) 


0.03764(23) 


2.30 


6 


-0.429(17) 


0.9526(11) 


0.03609(37) 


0.51 


8 


-0.412(47) 


0.9520(19) 


0.03591(58) 


0.54 


10 


-0.32(11) 


0.9501(28) 


0.03535(84) 


0.49 


12 


-0.42(17) 


0.9517(34) 


0.03581(98) 


0.50 



Table 7: Results from fitting the susceptibilities at j3 = 0.383245 
to eq. (|). 



Lmin 


c 


d 


V 


X 2 /dof 


4 


-0.4631(95) 


0.9502(10) 


0.03528(35) 


0.27 


6 


-0.551(93) 


0.9520(24) 


0.03580(72) 


0.22 


8 


-0.56(12) 


0.9522(26) 


0.03587(77) 


0.24 


10 


-0.58(15) 


0.9524(31) 


0.03594(93) 


0.26 


12 


-0.57(14) 


0.9523(29) 


0.03589(86) 


0.28 



Table 8: Fitting the susceptibility at fixed Q = 0.6240. 



Our fit results are given in table § and plotted in figure ||. As a final estimate 
for rj we quote 

^ = 0.0359(10). (13) 



Conclusions 

We have demonstrated that the spin-1 Ising model with suitably chosen cou- 
pling constants has remarkably improved finite size scaling properties. This 
allowed us to obtain high precision estimates for the critical indices v and r\ 
and two other universal quantities, the Binder cumulant Q and the ratio of 
partition functions Z a /Z p . 

Our results compare very well with a selection of estimates obtained by 
other authors, quoted in table || Many more estimates can be found in 
refs. H and f5|. Our estimate for the cumulant, Q* = 0.6240(2) may be 
compared with the one given in ||, namely Q* = 0.6233(4). 

It would be worthwhile to use the present model in studies of physical 
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0.039 I 1 1 1 1 1 1 

0.038 - 
0.037 

0.036 r. I T : 

0.035 - I 

0.034 - 

0.033 ' 1 1 1 1 1 

4 6 8 10 12 

^min 

0.039 I 1 1 1 1 1 1 

0.038 - * r 

0.037 I I | ] 

j. i * ] 

0.036 -. j 

0.035 - f ......... .. I - 

0.034 - 

0.033 ' 1 1 1 1 1 

4 6 8 10 12 

^min 

Figure 5: r] estimates from fit function eq. (f|) as function of L min . 
Top: at our estimated value for (3 C = 0.383245. Bottom: at (3 cho- 
sen two standard deviations above (upper data) and below (lower 
data) our f3 c estimate, i.e. at /3 = 0.383243 and f3 = 0.383247. 
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0.039 | 1 1 1 1 r- 

0.038 - 

0.037 - 
0.036 

0.03? I 1 1 • 

0.034 - 

0.033 I 1 1 1 1 ' 

4 6 8 10 12 

^min 



Figure 6: rj estimates as function of L min , fitting %-data at fixed 
Q to eq. (0). 



Ref. 


Method 


V 


rj 


i 


3D FT 


0.6301(5) 


0.0355(9) 


i 


3D FT 


0.6304(13) 


0.0335(25) 


§ 


EPS 


0.6293(26) 


0.036(6) 


E 


HT 


0.6310(5) 




E 


MC 


0.6301(8) 


0.037(3) 


i 


MC 


0.6308(10) 




present work 


MC 


0.6299(3) 


0.0359(10) 



Table 9: Comparing the results of the present study with pre- 
vious estimates in the literature. 3D FT means field theoretic 
calculations directly in three dimensions, EPS refers to e expan- 
sion, and MC means that the estimate relies on Monte Carlo 
simulations. 
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quantities not discussed in this work and to check to what extent the im- 
proved scaling behaviour helps to get better estimates. 

Of course, one could also search for even further reduction of corrections 
to scaling by refining our procedure or choosing variants of the model, e.g., 
the 4 -model or an Ising model with more couplings. Last but not least, 
application of the ideas underlying the present analysis to other models seems 
promising. 
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